Quantum Phases of Attractive Matter Waves in a Toroidal Trap 
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Investigating the quantum phase transition in a ring from a uniform attractive Bose-Einstein con- 
densate to a localized bright soliton we find that the soliton undergoes transverse collapse at a critical 
interaction strength, which depends on the ring dimensions. In addition, we predict the existence 
of other soliton configurations with many peaks, showing that they have a limited stability domain. 
Finally, we show that the phase diagram displays several new features when the toroidal trap is set 
in rotation. 



PACS Numbers: 03.75.Kk 
I. INTRODUCTION 

Metastable states of Bose Einstein condensates (BECs) 
made of attractive ^Li atoms, namely single and multi- 
ple bright soliton configurations, have been observed in 
two different experiments [1,2]. These metastable bright 
solitons, which can travel for long distances without dis- 
persion, have been the subject of various theoretical in- 
vestigations because of their relevance in nonlinear atom 
optics [3,4]. Recently, repulsive BECs in a quasi-lD ring 
have been produced and studied [5]. The case of an at- 
tractive BEC in a ring has not yet been experimentally in- 
vestigated but appears very interesting. For this system 
a quantum phase transition from a uniform condensate to 
a bright soliton. has been predicted by G.M. Kavoulakis 
and by R. Kanamoto, H. Saito and M. Ueda [6]. This 
prediction is based on mean-field and beyond mean-field 
numerical results for a ID Bose gas with contact inter- 
action and periodic boundary conditions [6]. Later, the 
same authors have shown that the quantum transition 
properties of the attractive BEC in a ID ring are strongly 
modified if the confining trap is rotating [7] . 

In this paper we investigate an attractive BEC in a 
3D ring, taking into account transverse variations of the 
BEC width, showing that the phase diagram of the sys- 
tem reveals novel and peculiar structures. In particular, 
we prove that, contrary to the simple ID case, the lo- 
calized soliton has a limited existence and stability do- 
main, which nevertheless strongly extends the stability 
domain of the uniform solution. Moreover, we find that 
the system supports also multi-peak solitons, which are 
energetically unstable but can be dynamically stable. Fi- 
nally, we analyze the effect of a rotating ring. In this case 
the multi-peak solitons are always energetically and dy- 
namically unstable, while the one-peak soliton is stable 
in a domain that, for a fixed rotation frequency, critically 
depends on the system parameters. 



II. TOROIDAL TRAP 

We consider a BEC confined in a toroidal potential 
given by 



(1) 



where p is the cylindric radial coordinate, ( is the cylin- 
dric axial coordinate and 9 is the azimuthal angle. The 
BEC has transverse harmonic confinement of frequency 
and the two characteristic lengths of the toroidal trap 
are Rq and a±_ — {h/{muj±))^^'^ . In the remaining part 
of the paper we use scaled units: time in units of 
length in units of a± and energy in units of huj±. To 
simplify the 3D Gross-Pitaevskii equation (GPE) we im- 
pose that the order parameter ^{p,9,C,,t) of the Bose 
condensate is the product of a generic azimuthal func- 
tion f{0,t) and a radial Gaussian function, which has 
two variational parameters: the width a{9, t) and the co- 
ordinate R of the center of mass in the radial direction. 
The trial wave function is given by 



f{9,t) ^^^y- 2.iitr 

i?l/2 TT^/^a{9,t) 



(2) 



We follow a procedure similar to that described in [8] in- 
serting the trial wave function into the GPE Lagrangian 
density 



4 + iv^-f/-imi^' 



(3) 



where F = AnasN/a^, with N the number of condensed 
atoms and Os < the attractive s-wave scattering length. 
The trapping potential of Eq. (1) has a cusp at the ori- 
gin, a feature not usually present in experimental traps, 
but this is expected to have minor effects on physical 
properties if only a small fraction of particles are near 
the origin. For example, if R/a > 3/2, the fraction of 
particles belonging to the radial region [0,R — a] is below 
1%. In the range 1 < R/a < 3/2 the population near the 
origin is not so small and the results reported below are 
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not fully reliable in this range. We integrate over p and 
C, coordinates [9]. In this way, from the Euler-Lagrange 
equations of the resulting effective Lagrangian density, 
we get the nonpolynomial Schrodinger equation (NPSE) 



.d_ 1 9^ 

' dt ^ 2 9^2 



-T{n) 



(4) 



where il){z^ 
tion of the 
sity profile 
W{n) = n{ 
is given by 
ometry the 
conditions 
variational 
is given by 



t) = f{e,t)/R^/'^ is the azimuthal wave func- 
condensate with z ~ R9, n — \ip\'^ is the den- 
normalized to unity and T{n) = '^'^^"•^ with 
1 — 5n)^/2. The scaled interaction strength g 
g = |r|/(27r) = 2N\as\/a^. In toroidal ge- 
solution t) must obey periodic boundary 
^(0,i) = ip{L,t) with L = 2itR. Within our 
approach the transverse width a of the BEC 



2 (^ \l/2 
(7 ^[l-gn) ' 



(5) 



For (/n <C 1 one has cr ~ 1, r(n) ~ —gn + 1, and the 
NPSE reduces to the ID GPE. In addition, we find that 
the variational parameter R is implicitly given by the 
equation 



R 



I 

R JO 



2-kR 



dz 



dz 



9 
2a2 



\il^\^dz = Ro . (6) 



This formula shows that the effective radius R of the 
BEC ring depends on the interaction strength g. We have 
verified that for a static and attractive BEC, the effective 
radius R decreases very slowly by increasing g, while for 
a repulsive BEC the opposite is true. In practice, because 
for an attractive BEC the strength g cannot exceed the 
scaled inverse density , the effective radius R is always 
close to Rq: for a uniform BEC, where Eq. (6) becomes 
R + g/ (47ri?2(l - g/{2nR)y/^) = Rq, it is easy to check 
that the relative difference between R and Rq is tipically 
of few percents and becomes 10% only near the collapse. 

The very good accuracy of the NPSE in approximat- 
ing the 3D GPE with a transverse harmonic potential has 
been verified in [8] for both positive and negative scat- 
tering length. In the derivation of NPSE one neglects 
the space and time derivatives of the width a{z,t). By 
including these terms one gets the coupled equations de- 
rived by Kamchatnov and Shchesnovich [10], but it is 
not clear if these terms give an improvement. We have 
verified that, according to the 3D GPE, the single-peak 
bright soliton of an attractive BEC in an infinite cylin- 
der collapses at the critical strength gc/2 = 0.676 (see 
also [11]), a value very close to the NPSE prediction: 
gc/2 = 2/3 = 0.666 [8]. 



III. UNIFORM AND LOCALIZED SOLUTIONS 

The NPSE conserves both the norm of the wavefunc- 
tion and the total energy E of the configuration. The 



stationary solutions follow from Eq. (4) by looking for so- 
lutions of the form '(p{z,t) = (/)(z)e~*^* for some chemical 
potential /j,. The resulting non-linear eigenvalue equation 
is the static NPSE. In toroidal geometry, the uniform so- 
lution (l){z) = I/Vl is always present for g < L, i.e. with 
density N/L < aj^/(2|as|), and corresponds to the eigen- 
value /i = r(l/L). In addition other, less trivial, profiles 
may be present. Beyond this limit (i.e. for g > L) the 
attraction is too strong and no regular solution is possi- 
ble, leaving the BEC collapse as the only possibility. A 
generic (real) solution (l){z) of the stationary NPSE may 
be interpreted as the classical "time" evolution of a ficti- 
tious particle moving in a potential V{(j)) = ^icfP — W {(j)^) . 
As a consequence, the "energy" conservation equation for 
this motion reads: 



1 /d4) 
2 



y(0) 



(7) 



According to the values of the two parameters /x and e two 
kinds of "trajectories" may be realized: For /x > and 
e > the solution (j){z) oscillates between a positive and 
a negative value, thereby crossing zero, while for e < 
the solution (j}{z) remains always positive. In the first 
case the solutions are named nodal solitons while in the 
second case nodeless solitons. 
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FIG. 1 (color online). Existence diagram in the {R,g) 
plane, where g — 2N\as\/a± is the interaction coupling 
and R = L/{2'k) is the azimuthal radius of the ring (in 
units a±). Uniform solution: exists for g < 2ttR. One- 
peak localized nodeless solution: exists between the two 
solid lines. Nodeless two-peak localized solution exists 
between the two dashed lines. The nodal two-peak local- 
ized solution exists between the two dot-dashed lines. 

For fixed g and L, the two parameters fi and e are 
implicitly determined by the two consistency equations: 



L = V2Ns 



1 



(8) 



and 
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rtpmaa 

V2N, / 

J'l'm.in 



d(t) 



(9) 



where Ng is the number of maxima of the soUton, i.e. 
the number of "periods" of the corresponding orbit, and 
(f'min i4>max) IS the minimum (maximum) value attained 
by <j){z) during the periodic oscillation. As customary in 
Newtonian problems, the extremal values of cj) are im- 
plicitly defined by the equation V{4>) = e in terms of e 
and n- The first equation comes from the commensura- 
bilty requirement of the solitonic solution, which, after 
an integer number of periods must close without discon- 
tinuities. The second equation is just the normalization 
condition of the wavefunction. 

The two consistency equations (8) and (9) have been 
solved and the domain of existence of soliton solutions of 
given topology has been determined in the {R, g) plane. 
Results are shown in Fig. 1. The numerical analy- 
sis shows that two solutions of same symmetry may be 
present in a portion of the existence domain. For sake of 
clarity in Fig. 1 we report the existence domain only for 
the one-peak soliton {Ng = 1) and for the two-peak soli- 
ton {Ng = 2). Note that the existence diagram of the ID 
GPE is strongly different from that of the NPSE. In fact, 
as previously stressed, the ID GPE does not take into 
account the transverse dinamics and, as a consequence, 
no collapse of solitonic solutions is predicted by ID GPE. 



IV. ENERGETIC STABILITY 

We now turn to the discussion of the energetic stabil- 
ity [12] of the previously defined solitonic configurations. 

One can rigorously prove that the energetic stability can 
be expressed in terms of the eigenvalues A; oi H ± nT', 
where 



H 



1 , dTin) 



and n(z) is the density profile of the stationary solution. 
The stationary solution (j){z) is energetically stable only 
if either of the two conditions is satisfied: 1) all the eigen- 
values A; are non-negative; 2) a single negative eigenvalue 
Ao is present and < 0. When we apply this general 
analysis to the simple case of the uniform stationary so- 
lution in toroidal geometry 0(2;) = we find that 
this case satisfies the latter of the two conditions previ- 
ously stated and the solution is stable until the second 
eigenvalue gets negative, triggering the instability. The 
resulting stability condition (27r/L)^-|-4nT' > explicitly 
becomes: 



5(' 



3/2 



> 1 



35 

4L 



(11) 



This formula reduces to -k"^ /{gL) for large L (ID limit), 
that is precisely the result one finds with the ID GPE 
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[6]. The stability analysis of the solitonic configurations, 
however, does not allow for a general analytic solution 
and the eigenvalue equations of if ± nT' must be inves- 
tigated numerically. 
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FIG. 2. Energy E and chemical potential [i of the one- 
peak solitons as a fimction of the coupling g for L = 10. 
The energetically stable solutions are the solid ones, the 
unstable solutions are the dashed ones. 

The operators H ± r\T' may be numerically diagonal- 
ized by introducing a finite mesh in the interval (] < z < L 
and approximating the differential operator with the cor- 
responding finite difference operator. The two resulting 
equations then give rise to a matrix eigenvalue problem. 
The numerical results show that: 

• in the regions where the uniform solution satisfies 
the energetic stability condition ^ < no other 
solitonic wavefunction can be stabilized; 

• only the one peak soliton is energetically stable in 
a portion of the domain where it is defined; 

• when distinct one peak solutions exist for the same 
values of g and L, the soliton is stable only in the 
branch corresponding to the lowest energy. 

The latter remark is illustrated in the upper panel of 
Fig. 2 where the energy E of the one-peak solution is 
shown as a fimction of g for L = 10. It is also interesting 
to plot the chemical potential /i versus g in the stable and 
unstable branch. The case of L = 10 previously analyzed 
is shown in the lower panel of Fig. 2, where it appears 
that the onset of instability corresponds to an extremum 
of the coupling constant g as a function of the chemical 
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potential. In fact, this immediately follows from the pre- 
vious analysis which led to ^ < 0: if H + nT' admits 

a single negative eigenvalue, ^ = signals the onset of 
the instability. 




FIG. 3 (color online). Attractive BEG in a ring ro- 
tating with frequency fl. Energetic-stability diagram in 
the plane {R, g) . The uniform solution is the ground- 
state only below the solid line. The one-peak localized 
solution exists between the two dashed lines but it is the 
ground-state only between the solid line and the upper 
dashed line. Note that for almost all R — L/(27r) the 
solid curve is superimposed to the lowest dashed curve: 
only for 1 < i? < 1.5 the lowest dashed line is below the 
solid line. 

The energetic stability region of stationary solutions of 
the NPSE in a ring are shown in the top panel (il = 0) 
of Fig. 3. Below the solid curve the uniform solution is 
energetically stable. The one-peak soliton, which exists 
between the two dashed curves, is energetically stable in 
the domain limited by the solid and the upper dashed 
line. In the remaining regions of the phase diagrams no 
energetically stable stationary solution is present and the 
BEG is expected to collapse. Note that for large R = 
L/{2tt) the upper dashed line tends to g = 4/3, that is 
the formula one finds for the collapse of a bright soliton 
in a infinite cylinder (see L. Salasnich et al. 2002 in [3]). 
It is not difficult to show that with a large R the existence 
domain of a A^s-peak bright soliton is instead given by 



0<g<-N,. 



(12) 



See for instance the existence domain of the two-peak 
solitons shown in Fig. 1. As previously stressed upper 
bounds do not exist within the ID GPE approach: the 
collapse of single and multiple bright solitons is due to 
the transverse dynamics of the condensate. 

In the other panels of Fig. 3 it is shown the effect 
of a rotating toroidal trap on the attractive BEG. The 
analysis is developed by observing that one has to include 
the centrifugal operator 



(13) 



into Eq. (4), where Q is the rotation frequency (in units 
of the frequency ll}± of the harmonic transverse confine- 
ment) and Mg is the azimuthal angular momentum. As 
shown in Ref. [7] by using the ID GPE, the uniform state 
of the attractive BEG is superfluid, i.e. it exists a criti- 
cal frequency flc below which the uniform state remains 
stationary, and only above this critical frequency the uni- 
form state rotates. In general the stationary uniform so- 
lution is thus given by ip{z) = e^^'^^-'/^/v'i, where the 
integer number j is a function of fl and L, namely 



j(ri, L) = int 



47r2 



(14) 



with int[x] the maximum integer that does not exceed 
X. The localized soliton solution has a different behav- 
ior: its angular momentum is not quantized [7] and this 
means that the quantum phase transition from the uni- 
form to the localized state suppresses the superfluidity 
of the system [7]. Setting ^{z) — 0(z)e*"^^\ where both 
(j){z) and a{z) are real, from the stationary NPSE with 
the centrifugal operator of Eq. (13) one finds 



da L c 

where the constant c is given by the equation 

dz 



+2^ 



0(Z)2 



] 



(15) 



(16) 



In addition, the function 4>{z) is obtained from the two 
consistency equations (8) and (9) with V{(t>) now given 

by 



202 ' 



(17) 



where W{(j)) = - 50^)1/2. 

Following the previous analysis one finds that the en- 
ergetic stability condition reads 



r>£2 

47r2 



7r2 (1 



£)3/2 



(18) 



In the ID limit of large L the previous formula gives 
1-4 - L)j > gL/n'^, that is the results found 
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in Rcf. [7]. Fig. 3 shows the energetic stabihty diagram 
in the plane (R, g) for two non-zero values of the rotating 
frequency Q. Below the solid line the stability condition 
holds and the uniform state is energetically stable. The 
periodic structure for f2 7^ is a consequence of the peri- 
odic quantization of the angular momentum j(f2, L). In 
particular, for a fixed ft, the solid line touches the hori- 
zontal axis = for discrete values oiR = L/ (27r), which 
correspond to jumps in the quantum number j. 
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FIG. 4 Angular momentum Me = {Me) and azimuthal 
width S — (z^)^/^ of the one-peak soliton as a function of 
the coupling g for = 0.1 and two values of i?: i? = 2.2 
solid line, ii = 2.3 dashed line. 

The one-peak solitonic solution exists between the two 
dashed lines and it is energetically stable between the 
solid line and the upper dashed line. Interestingly, Fig. 

3 shows that for large R the lower dashed line and the 
solid line practically coincide. As in the non rotating 
case, we find that also for O ^ solutions with more 
than one peak are not energetically stable. 

In Fig. 4 we plot the angular momentum Mg and the 
azimuthal width E of the rotating one-peak bright soli- 
ton. The figure shows that the angular momentum is not 
quantized: it approaches the value = QR^ of a "clas- 
sical particle" for g close to the collapse {g ~ 4/3), but it 
becomes quantized, i.e. Mg = j, where j depends on Q 
and R, for the small value of g that gives the transition 
to the uniform solution. Obviously, when the angular 
momentum becomes quantized the width S of the bright 
soliton coincides with that of the uniform solution. Fig. 

4 also shows that the width E is independent on the ring 
radius R as the bright soliton is close to the collapse; in 
this case g ~ 4/3 and S ~ 0.85. 



V. DYNAMICAL STABILITY 

It is important to stress that energetic stability implies 

dynamical stability but the converse is not true [12]. In 
order to investigate the dynamical stability of station- 
ary solutions in the ring one can solve the Bogoliubov-de 
Gennes (BdG) equations which give the elementary exci- 
tations ei {I = 1, 2, 3, ...) of the system: the appearence of 
a complex excitation signals dynamical instability while a 
negative excitation implies energetic instability [12]. For 
the uniform solution tjj{z) = e'^'^^^/^/VX one finds 



47r2 

n-j^j{n,L)]i 




4.9 (i + lf) 

L (l-f)3/2 



1/2 



(19) 



This result confirms that the uniform solution is energeti- 
cally stable if the previously written stability condition is 
satisfied; moreover it shows that the dynamical stability 
of the uniform solution is independent on fl. For localized 
solutions the BdG equations are computationally rather 
demanding. For this reason we have analyzed the dynam- 
ical stability by numerically solving the time-dependent 
NPSE taking as initial condition the stationary localized 
solution (j){z) with a very weak perturbation. 
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FIG. 5. Left panels: density profile p{z) of the solitonic 
solutions. Right panels: time-dependence of the mean 
squared widths (z^) — (z)^ for the weakly perturbed soli- 
tonic solutions. Ring axial Icnght: L = 15. Interaction 
strength: g = I. Rotational frequency: f2 = 0. 

In Fig. 5 we plot the density profile p{z) = |(/>(2;)|^ of 
the one-peak and the nodal two-peak solutions, choosing 
L = 15, g = 1 and two values of Cl. In addition, we plot 
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the time-evolution of the mean squared width (z^) — (z)^. 
Its behavior reveals that these solitonic sohitions are dy- 
namically stable. We have investigated the dynamical 
stability for various initial conditions. For f2 = we 
have found that: 

• the one-peak soliton is dynamically stable where it 
exists; 

• the nodal two-peak soliton is dynamically stable 
in the plane {R, g) only below the upper curve of 

existence of the one-peak soliton; 

• the nodeless two-peak soliton is dynamically unsta- 
ble. 

Similar results are found for solitonic solutions with a 

larger number Ns of peaks. 

The case ^ leads to similar results, keeping in 
mind, however, that nodal solitons do not exist under ro- 
tation. For high rotational frequencies, namely when f2 
approaches the frequency of transverse harmonic confine- 
ment (that is 1 in our units), the effect of the centrifugal 
force on the transverse dynamics becomes relevant. For a 
given angular momentum j ~ iJ^Q, the centrifugal force 
increases the cifFcctivc; radius R of the BEC ring. For 
a rotating and uniform ideal BEC one easily finds from 
Eq. (6) that R = Ro/{l - fi^). This formula holds also 
for an azimuthally localized ideal BEC. At the critical 
frequency = 1 the radius R diverges and this means 
that the Bose condensate is no more confined. As in the 
non-rotating case, an investigation of Eq. (6) shows the 
effect of the interaction strength g on the effective radius 
R is negligible for an attractive BEC. 

VI. CONCLUSIONS 

We have predicted novel quantum phases for an attrac- 
tive Bose condensate in a ring. Our results can be exper- 
imentally tested with the optical and magnetic traps re- 
cently developed [5,13], where a time-dependent stirring 
potential can be used to set into rotation the system. For 
instance, by choosing 10^ ^Li atoms (a^ = — 1.4 nm) in a 
toroidal trap with L ~ 25 /xm and ox — 3 /um, a sequence 
of transitions between the uniform state and solitonic 
configurations takes place when the stirring frequency 
CI is varied between and 1 kHz. These experiments 
will open the way to the observation of amazing quan- 
tum phenomena like the solitonic condensate without su- 
perfluidity and the dynamically induced phase transition 
from uniform to localized states. 
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